function ecog = ecogBandpassFilter(ecog, band, banddefinition)

% this function filters the data of each condition in a desired frequency band.
%
% INPUT:  (1)  the ecog structure containing the fields "conddata"
%         (2a) a vector of numbers between 1 & 6 together with 'predef'
%         (2b) the lower and upper limit of the band pass filter together
%              'selfdef'
%         (3)  banddefinition: either 'predef' for the predefined bands in
%              ecog.freqband or 'selfdef' for a selfdefined band within an
%              arbitrary range.
%
% OUTPUT: (1)  ecog structure with the new fields
%              "freqband" ==> each row contains the lower and upper limit
%                             of the frequency band
%              "bandname"
%              "filtdata" ==> filtered trials of each condition and each
%                             band
%              "filtavr"  ==> filtered average of each condition and band
%
% EXAMPLE: ecog = ecogBandpassFilter(ecog, 1:6, 'predef') for all predefined 
%          frequency bands.
%          ecog = ecogBandpassFilter(ecog, [5 40], 'selfdef') for a
%          selfdefined band between 5 and 40 Hz.
%
% REVISION: 06/25/2010; SD wrote it

clc;


if isfield(ecog,{'filtdata', 'filtavr'}) > 0; 
    ecog = rmfield(ecog,'filtdata');
    ecog = rmfield(ecog,'filtavr');
end


switch banddefinition
    case 'predef'
        if ecog.srate < 400;
            ecog.freqband = [4 8; 8 10; 10 13; 13 30; 30 80; 80 110]; 
        else
            ecog.freqband = [4 8; 8 10; 10 13; 13 30; 30 80; 80 200]; 
        end
        ecog.bandname = {('theta'),('lowalpha'),('uppalpha'),('beta'),('gamma'),('higamma')};
        ecog.lowfreq  = ecog.freqband(band,1);
        ecog.uppfreq  = ecog.freqband(band,2);
        
    case 'selfdef'
        ecog.freqband = band;
        ecog.lowfreq  = ecog.freqband(1);
        ecog.uppfreq  = ecog.freqband(2);
                
end


for k = 1:size(ecog.freqband,1);                                                          
    if size(ecog.data,3) == 1;
        [b,a] = butter(4, [ecog.lowfreq(k)/(ecog.srate/2), ecog.uppfreq(k)/(ecog.srate/2)]);
        ecog.BPdata{k} = filter( b,a,ecog.data,[],2);
    elseif size(ecog.data,3) > 1;
        % this is the case after artifact rejection
        [b,a] = butter(2, [ecog.lowfreq(k)/(ecog.srate/2), ecog.uppfreq(k)/(ecog.srate/2)]);   
        for trial = 1:size(ecog.data,3);
            ecog.BPdata{k}(:,:,trial) = filter( b,a,ecog.data(:,:,trial),[],2);
        end
    end
end




if isfield(ecog,'conddata') == 1;
    for cond = 1:size(ecog.conddata,1);
    for k = 1:size(ecog.freqband,1);                                                
        [b,a] = butter(4, [ecog.lowfreq(k)/(ecog.srate/2), ecog.uppfreq(k)/(ecog.srate/2)]);

        Trials1 = zeros(size(ecog.conddata{1,2}));                      
        Trials2 = zeros(size(ecog.conddata{1,1}));

        for trial = 1:size( ecog.conddata{cond,2}, 3 ); 
            Trials1( :,:,trial ) = filter( b,a,ecog.conddata{cond,2}( :,:,trial ),[],2 );
        end
        ecog.filtdata{cond, k} = Trials1;
        for trial = 1:size( ecog.conddata{1,1}, 3 ); 
            Trials2(:,:,trial) = filter( b,a,ecog.conddata{cond,1}( :,:,trial ),[],2 );
        end
        ecog.filtavr{cond, k} = Trials2;
%         fprintf('cond: %i -- band number: %i\n', cond, k)
    end
    end
end
end